Recent advances in nanomagnetism technology allow the fabrication of Stoner particles (single-domain magnetic nanoparticles due to strong exchange interaction),^[1–4] which has attracted much attention both from the fundamental point of view and from the potential applications in information industry.^[5,6] Magnetization dynamics of a single Stoner particle has been extensively studied, either by static or time-dependent magnetic field method^[7–17] or by spin-polarized current approaches.^[18–23] The current-induced magnetization switching attracts much attention due to the locality and convenient controllability of current. However, the high critical switching current I_c (or density) hinders practical applications.^[24] Many efforts were devoted to lowering the current I_c, for instance, by designing an optimized spin current pattern,^[23] by using pure spin current^[25] and thermal activation.^[26,27] Recently lower critical current has also been proposed in ferromagnetic semiconductors.^[28–30]

Moreover, since magnetic nanoparticles were actually fabricated into array,^[1–4] the dipole–dipole interaction (DDI) between them becomes important for the magnetic switching behavior. Indeed there were several studies on systems of two Stoner particles, the simplest case where the DDI effect is investigated, concluding that the dipolar interaction can assist the magnetic switching.^[31–37] Recently we revisited the problem of two-body Stoner particles and proposed a novel technological perspective by showing that in a synchronized dynamic mode (i.e., both the magnetic vectors of the particles run synchronously), the critical switching field in presence of magnetic DDI can be dramatically lowered. This implies that even zero-field switching can be achieved by appropriately engineering the DDI strength through adjusting the interparticle distance.^[38] In the present study, we extend these studies by investigating the current-induced synchronized magnetization switching of such a two-body Stoner particle system on injecting a spin-polarized current through a spin-valve-like structure. As a numerical result, the critical switching current I_c (or its density) can be lowered to about 70% of the usual value without DDI around a critical DDI strength. As is well known, high switching current is a major hindrance to the direct application of current-induced magnetization switching, our results may have an important reference for its future industrial applications. It is noticed that the critical DDI strength corresponds to a distinguished interparticle distance, which is the same value as that in the field-induced magnetization switching.^[38] The critical current I_c bifurcates at a critical DDI strength by a square-root-like behavior. Moreover, we numerically investigate the issues including the magnetic hysteresis loop, magnetization self-precession, reversal time, and synchronization stability phase diagram for the two-body Stoner particle system in the presence of the DDI. These investigations reveal many novel promising phenomena relating to future device applications. The rest of this paper is organized as follows. In Section 2, the two typical geometrical configurations of two-body system in a spin-valve-like setup are introduced along with the equations governing current-induced magnetization dynamics. Our analytical and numerical results are presented in Sections 3 and 4, respectively. Finally, some conclusions are drawn from the present study and also some discussion is presented in Section 5.

We consider the spin-valve structure schematically shown in Fig. 1. For simplicity, we consider two particles that are in the same sphere shape. Two typical geometric configurations of the two-body system are investigated, where the unit vector n along the line connecting the two particles is either perpendicular or parallel to the magnetic easy axis (EA), which is assumed to be along z axis. These two configurations are referred to as PERP configuration and PARA configuration, respectively. That is, if we define n = (sin θ_n cos ϕ_n, sin θ_n sin ϕ_n, cos θ_n) in the usual spherical coordinates where z axis is the pole axis, thus θ_n = π/2 for PERP configuration and θ_n = 0, π for PARA (the two values are equivalent since DDI is invariant under n → −n). Moreover, one can let ϕ_n = 0 for both configurations, i.e., n = x in PERP and n = z in PARA (see Fig. 1).

Fig. 1.

	Figure Option ViewDownloadNew Window
	Fig. 1. Schematic diagram of two typical geometric configurations (PERP and PARA) in spin-valve structure and corresponding coordinate axis orientations. Polarizer is a thick pinned ferromagnetic layer where s denotes its polarized direction and NM denotes nonmagnetic metallic layer. EA denotes easy axis and e denotes electrons flow.

The magnetization dynamics of two-body Stoner particles under a spin-polarized current is governed by the modified Landau–Lifshitz–Gilbert (LLG) equation^[39] with the spin-transfer torque (STT) term,^[19–21] 1 Here m_i = M_i/M_s is the normalized magnetization vector of the i-th-particle (i = 1, 2), M_s = |M_i| is the magnetization saturation of either particle, and α is the phenomenological Gilbert damping coefficient. Now, we assume that the two particles are completely identical, all with the same concerned parameters for simplicity. The time is set to be in 1/(|γ|M_s) where γ is the gyromagnetic ratio. The total effective field on each particle comes from the variational derivative of the total magnetic energy with respect to magnetization, . For the concerned system, the energy E per particle volume V (in units of , with μ₀ being the vacuum permeability) can be expressed as 2 Here the uniaxial magnetic anisotropic parameter k covers both shape and magnetocrystalline contributions to the magnetic anisotropy along the easy axis. The parameter η = V/4πd³ is a geometric factor characterizing the DDI strength, where d is the fixed distance between the two particles. The exchange interaction energy between two particles is omitted in Eq. (2) since it becomes important only at very small particle distances. Moreover, in the synchronized magnetic dynamics to be investigated below, it only contributes a constant to the energy and will therefore not change the physical behavior. The LLG equation can be written in the usual spherical coordinates, by using m_i = (sin θ_i cos ϕ_i, sin θ_i sin ϕ_i, cos θ_i).

Spin transfer torque^[19–21] is an effect of angular momentum transfer from the conduction electron spins to local magnetic moment. The transferred angular momentum will generate a torque on the magnetic moment, driving it to reverse or precession due to the principle of angular momentum conservation. An important STT term T_stt = −a_Im_i × (m_i × s) is used (i = 1,2), which takes the Slonczewski’s form^[19,20] where has a (normalized) magnetic field dimension. The field-like STT^[21] is omitted because it is usually much smaller than the Slonczewski’s term. Here (among standard notation) I and P are the magnitude and spin polarization degree of the current, respectively, s is the unit vector along the current polarization direction (see Fig. 1) which is characterized by two angles θ_s and ϕ_s in the spherical coordinates. Since a_I is proportional to the current I (assuming other parameters to be constant), we will refer to a_I as the switching current in the following description. For instance, a_I,c denotes the critical switching current.

In the following, we will concentrate on the synchronized motion mode of the two dipoles in the current-induced switching. Both magnetization vectors remain parallel throughout the motion, θ₁ = θ₂ = θ, and ϕ₁ = ϕ₂ = ϕ. Analytical results can be obtained in such a motion and the stability of the synchronization will be numerically verified later. The two-body Stoner particles promise to act as an information bit in computer technology. The modified LLG equations (1) in spherical coordinates are a system of coupled nonlinear equations and read 3 4 Here ψ_n/s is the angle between m and n/s such that cos ψ_n/s = cos θ cos θ_n/s + sin θ sin θ_n/s cos (ϕ − ϕ_n/s). Note that we have defined ϕ_n = 0).

In order to find the critical switching current a_I,c, one should find the fixed points (FPs) of the nonlinear Eqs. (3) and (4), and then analyze their stability conditions.^[40] We first examine the FP condition for the two poles θ_n = 0, π since their stability analysis is directly related to the interested critical switching current a_I,c. From Eqs. (3) and (4), we obtain 5 where ∓ or ± corresponds to θ = 0 and θ = π, respectively. Due to the uncertainty of ϕ at poles, we conclude only sin 2 θ_n = 0 and sin θ_s = 0 such that the two poles become the possible FPs. Thus, we only need to focus on the simple cases of θ_s = 0, π and θ_n = 0, π, π/2 (i.e., PARA and PERP configurations). In the following, we let θ_s = π, which is the equivalent case for θ_s = 0 by inverting the current flow direction.

We first recover the usual case without DDI (i.e., η = 0), where PERP and PARA configurations are indistinguishable. The motion of θ is decoupled to ϕ, 6 where we define Γ = 1 + α². The critical switching current a_I,c can be obtained by analyzing the linear stability conditions at poles.^[40] One can find that the FP of θ = 0 becomes unstable but θ = π is stable when a_I > 2 αk. It implies that the magnetization reversal switches on and a_I,c reads 7 where the superscript index 0 denotes the zero DDI case. It is pointed out that 2k is the usual Stoner–Wohlfarth limit^[7] in the field reversal case, namely the barrier height between the two minima at poles in the energy landscape.

We now turn to the case with DDI (η ≠ 0). The magnetic energy (without external field) under the synchronized motion mode is 8 which leads to the equilibrium (minimal) energies for the PARA and PERP configurations as 9 Figure 2 shows that the essential differences in energy landscape between the PARA and PERP configurations. Figures 2(a)–2(c) show the energy variations in PARA configuration with the increase of the DDI strength η. The two edges of θ = 0, π should be understood as the two poles in a sphere. In the PARA case, the energy diagrams have no essential variations of geometrical shape with η. The two poles (θ = 0, π) are always the energy minima, and the in-between barrier (at θ = π/2) height increases with η, expecting a higher critical switching current (or field) in the PARA case.

Fig. 2.

	Figure Option ViewDownloadNew Window
	Fig. 2. (color online) Three-dimensional energy diagrams for the PARA and PERP configurations. ((a)–(c)): PARA with η = 0.1, 1/3, 0.9; ((d)–(f)): PERP with η = 0.1, 1/3, 0.9. System parameter k = 0.5.

However, in the PERP configurations as shown in Figs. 2(d)–2(f), the energy landscapes exhibit an essential geometrical change with the increase of η. Remarkably, it follows from Eq. (9) that a critical DDI strength η_c = 2k/3 exists, which is the value shared both in the current case and in the field switching case. Furthermore, η_c is equivalent to the interparticle distance where is the standard anisotropy coefficient.^[38] At the critical distance, the energy diagram has two flat paths along ϕ = 0, ϕ [see Fig. 2(e)], which means that any position of θ can be in equilibrium. While on the both sides of the critical η_c value, the system synchronized ground states (i.e., energy minimum) change. When the interparticle distance d > d_c (η < η_c), the ground states are at two poles of θ = 0, π [Fig. 2(d)] while when d < d_c(η > η_c), the ground states are at other two positions, θ = π/2 and ϕ = 0, π, along the hard x axis.

For the PARA configuration, one can find the decoupled equation of θ to ϕ, 10 With the linear stability analysis, one can derive that the FP of θ = 0 becomes unstable but θ = π is still stable when a_I > α(2k + 3η). Hence, the critical switching current in the PARA case reads, 11 Thus, is always higher than for η ≠ 0, as it is expected that the barrier height between two energy minima always increases with η increasing.

The PERP configuration is much more interesting and exhibits novel physical behavior due to the geometrical structure change in the energy diagrams as shown in Figs. 2(d)–2(f), in which a novel zero-field switching mechanism was proposed.^[38] In PERP, the nonlinear equations (3) and (4) can no longer be decoupled between θ and ϕ, 12 13 Again according to the linear stability analysis approach, the linearized matrices at the FPs of θ = 0, π can be found when with 14 15 Here, diag denotes the diagonal matrix and the superscript (subscript) sign refers to θ = 0(θ = π). One can find that for a_I > α(2k − 3η/2), θ = 0 becomes unstable but θ = π is still stable, such that the critical switching current reads, 16 The obtained is applicable only for the regime of η < η_c because the two-body system in the PERP case changes its energy equilibrium for η > η_c, where θ = π/2 may become the FP as shown in energy diagram Figs. 2(d)–2(f). Thus, one requires to analyze the stability condition around the new FP for the η > η_c regime, which can be obtained under the current case to be 17 Then, the linearized matrix at this FP is 18 One can find that the FP becomes unstable but the pole θ = π is still stable when . Thus the critical switching current within the η > η_c regime reads 19 So far we have obtained the analytical solutions of the critical switching current in both the PARA and PERP configurations, , , and , which are plotted by the solid, dashed and dot-dashed curves in Fig. 3, respectively, with the system parameters k = 0.5 and the damping α = 0.1. In principle, the minimum switching current in PERP should occur at the critical DDI value η_c, 20 which is obtained by using Eq. (16) at η_c and is half of the value without DDI. However, as we will show later [see Fig. 5(a)], there is a very small region around η_c above the crossed two solutions Eqs. (16) and (19), where the magnetization dynamics is not static but shows self-precession. Thus, in practice, we numerically find that the minimal critical switching current for a static magnetization reversal can be lowered by about 70% of the usual value without DDI, which will be more discussed later in Fig. 5(a).

Fig. 3.

	Figure Option ViewDownloadNew Window
	Fig. 3. (color online) Normalized critical switching current aI,c versus the DDI strength η for the PERP and PARA configurations. Curves show analytical results, and circles and squares denote numerical results. System parameters are k = 0.5 and α = 0.1.

In order to complement and quantitatively support our previous discussion, we now numerically verify our analytical results and investigate the reversal phenomena of two-body Stoner particle system in the presence of magnetic dipolar interaction under a spin-polarized current. We solve LLG Eq. (1) by implementing a fourth-order Runge–Kutta method. We consider the range of the DDI strength to be 0 ⩽ η ⩽ 1. In this case, η = 0 corresponds to the limiting case of two particles infinitely apart from one another i.e., d → ∞. Large η may be realized by fabricating magnetic nanoparticles of ellipsoidal shapes allowing closer proximity to one another. Throughout the numerical calculations, we use the damping parameter of α = 0.1. In Fig. 3, we first compare the numerical simulations of the critical switching current a_I,c between PERP and PARA configurations (circles and squares) and show excellent agreement with our analytical solutions: Eq. (11), Eq. (16), and Eq. (19). The critical switching current in PARA is always higher than the value without DDI (at η = 0), while in PERP it bifurcates at the critical DDI value (at η = 1/3) with a square-root dependence behavior.

4.1. Hysteresis loop

With the numerical simulations, we study the magnetic hysteresis loops under the current-induced synchronized switching. Figure 4(a) shows the hysteresis loops at different DDI strengths for the PARA configuration. The increase of η will enlarge the loop size, implying large critical switching current. Figure 4(b) shows the hysteresis loops for the PERP configuration for η < η_c. The increase of η will reduce the loop size, in contrast with the PARA case. In particular, when η approaches to η_c, magnetic self-precession phenomenon in current switching can be observed, which will be discussed in detail in the next subsection. For example, there are very narrow regimes at η = 0.31, which is enlarged in the inset of Fig. 4(b), showing the m_z oscillation amplitude. Figure 4(c) shows the hysteresis behavior for the PERP configuration for η > η_c. It is interesting to note that there are no loops but an invertible tristate occurring in Fig. 4(c). Finally, for comparison, we present Fig. 4(d) to show the hysteresis behaviors in the field switching case for the PERP configuration, where h is the normalized external field by the unit magnetization saturation M_s. When η < η_c a clear hysteresis loop occurs while when η > η_c no hysteresis loops are found but the magnetization m_z is linearly stable with respect to h, different from the tristate phenomenon in the current driven case.

Fig. 4.

	Figure Option ViewDownloadNew Window
	Fig. 4. (color online) Hysteresis loops. (a) PARA configuration. (b) PERP configuration when η < ηc = 1/3. Insets show magnified magnetic self-precession regimes, where mz is proportional to oscillation amplitude. (c) PERP configuration when η > ηc = 1/3. (d) PERP configuration in the field switching case. Other parameters are the same as those in Fig. 3.

Fig. 5.

Figure Option ViewDownloadNew Window

Fig. 5. (color online) (a) Phase diagram in the PERP configuration of the magnetization self-precession in parameter η–aI plane, where the color is proportional to the self-oscillation amplitude and the two dashed lines show the analytical solutions of the critical switching current. (b) Self-oscillation curve of mz (t) at η = 0.31 and aI = 0.06. Inset shows the spatial trajectories of the magnetization vector and its projections. Other parameters are the same as those in Fig. 3.

4.3. Reversal time

The numerical results for the reversal time T_r versus switching spin-polarized current magnitude a_I at different DDI strengths for both PARA and PERP configurations are shown in Figs. 6(a) and 6(b). The reversal time is numerically determined from m_z = + 1 to m_z = −1 for either particle, which depends on the initial direction of the magnetization vector due to the metastability of the pole.^[38] Thus, in practice, we deviate the initial angles θ₁ (0) = θ₂ (0) = 0.001 in order to find a finite reversal time. The same initial values are of importance for the stability of the synchronized motion mode, which will be discussed in next subsection. Figures 6(a) and 6(b) show that in both PARA and PERP cases, the curves of reversal time versus current have similar variation trends, approximately inversely proportional to the current magnitude a_I for large currents. In the PERP configuration, the T_r behavior is different from that in the field switching case when the DDI parameter approaches to its critical value. A substantially shorter reversal time can be found to be around zero field regime in the field case^[38] when η ≈ η_c. In contrast, the reversal time curve is non-smooth in the current case when η is around its critical value as shown in the inset of Fig. 6(b).

Fig. 6.

	Figure Option ViewDownloadNew Window
	Fig. 6. (color online) Plots of reversal time Tr versus current magnitude aI. (a) PARA configuration: different DDI strength η = 0,0.1,0.5,1; (b) PERP configuration: η = 0,0.1,0.2,0.33,0.4,0.7. The inset shows magnified part of plot for η = 0.33 where non-smooth Tr is observed. Other parameters are the same as those in Fig. 3.

4.4. Stability phase diagram

As an important issue, we now examine the stability of the synchronized magnetization dynamic scenario in the PERP configuration based on Eq. (1). We will examine the stability from two aspects: the difference in initial angle between the two magnetization vectors and the difference in material parameter between the two particles, which may result in the instability and spoil the synchronization.

We first numerically examine the average stable magnetization m_1z + m_2z/2 (the value −1 means stable while 0 unstable) dependent on the deviation angle δ between the initial directions of the two magnetization vectors. In detail, the initial direction of one particle is fixed along the z axis, i.e., θ₁ = 0.001 and the other one is changed as θ₂ = 0.001 + δ where 0.001 is for eliminating the singularity in pole. We can plot the stability phase diagram in the a_I–η space as shown in Fig. 7(a) for the current case, and plot in Fig. 7(b) for the field switching case for comparison in the h–η space. The color difference is proportional to maximal δ_m (in units of degree) for guaranteeing the stability where the gray color denotes δ_m >1° in Fig. 7(a) and δ_m > 10° in Fig. 7(b). When δ > δ_m, the final average magnetization is found to be 0, implying that the synchronized dynamics is destroyed. There are two points from Figs. 7(a) and 7(b), which should be mentioned. (i) In both current and field switching cases, the important stability regimes are located just above the critical switching current or field lines shown in dashed or dot-dashed lines in corresponding figures, especially for η < η_c. “Islands” appear in such regions. (ii) It has a wider initial deviation angle limit δ_m for the synchronization stability in the field case (< 10°) than in the current case (< 1°). We also plot the case of η = 0.23 and a_I = 0.07 in Fig. 8, which is 70% of the usual switching current without DDI. Figure 8(a) shows the time evolutions of all the components of the magnetization vector. Figure 8(b) shows the maximal initial deviation angle for the synchronization stability. Consequently, we remark that the practical critical switching current for a static magnetization reversal can be lowered to about 70% of the usual value .

Fig. 7.

Figure Option ViewDownloadNew Window

Fig. 7. (color online) Synchronization stability phase diagrams in the PERP configuration for the initial direction deviation δ of the magnetization vectors in (a) current case and (b) field case. The color difference is proportional to maximal stability δm in units of degree where the gray color denotes δm > 1° in panel (a) and δm >10° in panel (b). Synchronization stability phase diagrams for different experienced (c) δaI’s in current case and (d) δh’s in field case. The color difference is proportional to the variation width or δhw where the gray color denotes and δ hw > 0.02. Other parameters are the same as those in Fig. 3.

Fig. 8.

	Figure Option ViewDownloadNew Window
	Fig. 8. (a) Plots of magnetization vector component versus time at η = 0.23 and aI = 0.07. (b) Plot of average stable magnetization m1z + m2z/2 versus initial angle deviation δ, where δm ≈ 1.2° is observed. (c) Plot of (m1z + m2z/2 versus current deviation δaI, where is found. Other parameters are the same as those in Fig. 3.

On the other hand, for two particles with slightly differences, for instance, in anisotropy K_i, volume V_i, and saturation magnetization (M_s,i (i = 1,2)), as discussed previously,^[38] the critical value of DDI becomes and the effective fields experienced by particles due to the external field or the STT here will be slightly different from each other. Hence, we also numerically examine the synchronization stability due to such a difference. In detail we fix the current torque on one particle to be a_1,I and change the other one as a_2,I = a_1,I + δa_I where δa_I denotes a small deviation for different currents. The field case is the same except introducing a small deviation δh = h₂ − h₁ for different external fields. We thus plot the stability phase diagram in the a_I–η space as shown in Fig. 7(c) for the current case, and Fig. 7(d) for the field case for comparison in the h–η space. The color difference is now proportional to the variation widths of and δh^w, in which the final average magnetization is found to be −1 for stability. The gray color denotes in Fig. 7(c) and in Fig. 7(d). From the figures, one can observe the stability “islands” structures above the critical switching current or field lines become smeary now, and the current switching case has a wider stability window due to the parameter differences between the two particles than that in the field case. Figure 8(c) shows an example for current difference on two particles in the case of η = 0.23 and a_I = 0.07, where the width for the average stable magnetization equaling −1 to sustain the synchronized motion mode.

In this work, we compare our results with a concrete magnetic material such as cobalt (Co) particles. The standard data are M_s = 1.4 × 10⁶ A/m, uniaxial strength K = 10⁵ J/m³, and α = 0.1.^[13] Thus such that the critical DDI strength η_c = 0.027. If we consider two spherical particles with the radius r, we have the DDI parameter η = r³/(3d³). Thus the critical DDI is reached at d_c = 2.3r. Assuming r = 10 nm and spin polarization P = 0.4, the critical switching “effective field” without DDI is Oe (1 Oe = 79.5775 A · m⁻¹), which translates the critical switching current density to about j₀ = 3 × 10⁸ A/cm². In our model with DDI we argue that the critical switching current can be lowered to 70% of the original value, i.e. about 2 × 10⁸ A/cm² when the two particles are engineered to be located near the critical distance d_c = 23 nm. Also in the case of Co, the time is in units of 1/(|γ|M_s) = 3.23 ps. From Fig. 6(b), the reversal time is infinitely long at the critical switching current point. Increasing the switching current will reduce the reversal time inversely. More importantly, the key issue for the technological aim in our two-body Stoner particle system is to maintain the synchronized motion mode against large deviation from initial conditions and/or material parameters. The stability region from the initial angle deviations in the current case (about 1°) is much smaller than that in the field case (about 10°). Thus the enhancement of the synchronization stability might be significant for future studies. Moreover whether there exists a zero switching current in contrast with the zero-field switching case with the aid of DDI is also an interesting issue.

In conclusion, we have investigated the magnetization reversal of two-body uniaxial Stoner particles in a stable synchronized motion mode, by injecting a spin-polarized current through a spin-valvelike structure. In the presence of magnetic dipolar interaction, the critical switching current for reversing the two dipoles is analytically obtained and numerically verified in two typical geometric PERP and PARA configurations. In the interesting PERP configuration, the critical switching current bifurcates at a critical DDI strength with a square-root behavior, near which it can be lowered to about 70% of the usual value without DDI. Moreover, we also numerically investigate the current-induced magnetization hysteresis loops, magnetic self-precession phenomenon, reversal time, and the synchronization stability phase diagram in the two-body system, which shows interesting predictions and is expected to be useful in future device applications.